The second edition of Think DSP is not for sale yet, but if you would like to support this project, you can buy me a coffee.
6. Discrete cosine transform#
The topic of this chapter is the Discrete Cosine Transform (DCT), which is used in MP3 and related formats for compressing music; JPEG and similar formats for images; and the MPEG family of formats for video.
DCT is similar in many ways to the Discrete Fourier Transform (DFT), which we have been using for spectral analysis. Once we learn how DCT works, it will be easier to explain DFT.
Here are the steps to get there:
We’ll start with the synthesis problem: given a set of frequency components and their amplitudes, how can we construct a wave?
Next we’ll rewrite the synthesis problem using NumPy arrays. This move is good for performance, and also provides insight for the next step.
We’ll look at the analysis problem: given a signal and a set of frequencies, how can we find the amplitude of each frequency component? We’ll start with a solution that is conceptually simple but slow.
Finally, we’ll use some principles from linear algebra to find a more efficient algorithm. If you already know linear algebra, that’s great, but I’ll explain what you need as we go.
Click here to run this notebook on Colab.
6.1. Synthesis#
Suppose I give you a list of amplitudes and a list of frequencies, and ask you to construct a signal that is the sum of these frequency components.
Using objects in the thinkdsp module, there is a simple way to perform this operation, which is called synthesis:
from thinkdsp import CosSignal, SumSignal
def synthesize1(amps, fs, ts):
components = [CosSignal(freq, amp)
for amp, freq in zip(amps, fs)]
signal = SumSignal(*components)
ys = signal.evaluate(ts)
return ys
amps is a list of amplitudes, fs is the list of frequencies, and ts is the sequence of times where the signal should be evaluated.
components is a list of CosSignal objects, one for each amplitude-frequency pair.
SumSignal represents the sum of these frequency components.
Finally, evaluate computes the value of the signal at each time in ts.
We can test this function like this:
from thinkdsp import Wave
amps = np.array([0.6, 0.25, 0.1, 0.05])
fs = [100, 200, 300, 400]
framerate = 11025
ts = np.linspace(0, 1, framerate, endpoint=False)
ys1 = synthesize1(amps, fs, ts)
wave1 = Wave(ys1, ts, framerate)
This example makes a signal that contains a fundamental frequency at 100 Hz and three harmonics (100 Hz is a sharp G2). It renders the signal for one second at 11,025 frames per second and puts the results into a Wave object.
Conceptually, synthesis is pretty simple. But in this form it doesn’t help much with analysis, which is the inverse problem: given the wave, how do we identify the frequency components and their amplitudes?
6.2. Synthesis with arrays#
Here’s another way to write synthesize:
def synthesize2(amps, fs, ts):
args = np.outer(ts, fs)
M = np.cos(PI2 * args)
ys = np.dot(M, amps)
return ys
This function looks very different, but it does the same thing. Let’s see how it works:
np.outercomputes the outer product oftsandfs. The result is an array with one row for each element oftsand one column for each element offs. Each element in the array is the product of a frequency and a time, \(f t\).We multiply
argsby \(2 \pi\) and applycos, so each element of the result is \(\cos (2 \pi f t)\). Since thetsrun down the columns, each column contains a cosine signal at a particular frequency, evaluated at a sequence of times.np.dotmultiplies each row ofMbyamps, element-wise, and then adds up the products. In terms of linear algebra, we are multiplying a matrix,M, by a vector,amps. In terms of signals, we are computing the weighted sum of frequency components.
The following figure shows the structure of this computation.
#TODO: Figure here
Each row of the matrix, M, corresponds to a time from 0.0 to 1.0 seconds; \(t_n\) is the time of the \(n\)th row.
Each column corresponds to a frequency from 100 to 400 Hz; \(f_k\) is the frequency of the \(k\)th column.
I labeled the \(n\)th row with the letters \(a\) through \(d\); as an example, the value of \(a\) is \(\cos [2 \pi (100) t_n]\).
The result of the dot product, ys, is a vector with one element for each row of M.
The \(n\)th element, labeled \(e\), is the sum of products:
And likewise with the other elements of ys.
So each element of ys is the sum of four frequency components, evaluated at a point in time, and multiplied by the corresponding amplitudes.
And that’s exactly what we wanted.
We can use the code from the previous section to check that the two versions of synthesize produce the same results.
ys2 = synthesize2(amps, fs, ts)
np.max(np.abs(ys1 - ys2))
np.float64(1.2789769243681803e-13)
The biggest difference between ys1 and ys2 is about 1e-13, which is what we expect due to floating-point errors.
Writing this computation in terms of linear algebra makes the code smaller and faster.
Linear algebra provides concise notation for operations on matrices and vectors.
For example, we could write synthesize like this:
where \(a\) is a vector of amplitudes, \(t\) is a vector of times, \(f\) is a vector of frequencies, and \(\otimes\) is the symbol for the outer product of two vectors.
6.3. Analysis#
Now we are ready to solve the analysis problem.
Suppose I give you a wave and tell you that it is the sum of cosines with a given set of frequencies.
How would you find the amplitude for each frequency component? In other words, given ys, ts and fs, can you recover amps?
In terms of linear algebra, the first step is the same as for synthesis: we compute \(M = \cos (2 \pi t \otimes f)\).
Then we want to find \(a\) so that \(y = M a\); in other words, we want to solve a linear system.
NumPy provides linalg.solve, which does exactly that.
Here’s what the code looks like:
def analyze1(ys, fs, ts):
args = np.outer(ts, fs)
M = np.cos(PI2 * args)
amps = np.linalg.solve(M, ys)
return amps
The first two lines use ts and fs to build the matrix, M.
Then np.linalg.solve computes amps.
But there’s a hitch. In general we can only solve a system of linear equations if the matrix is square; that is, the number of equations (rows) is the same as the number of unknowns (columns).
In this example, we have only 4 frequencies, but we evaluated the signal at 11,025 times. So we have many more equations than unknowns.
In general if ys contains more than 4 elements, it is unlikely that we can analyze it using only 4 frequencies.
But in this case, we know that the ys were actually generated by adding only 4 frequency components, so we can use any 4 values from the wave array to recover amps.
For simplicity, I’ll use the first 4 samples from the signal.
Using the values of ys, fs and ts from the previous section, we can run analyze1 like this:
n = len(fs)
amps2 = analyze1(ys1[:n], fs, ts[:n])
amps2
array([0.6 , 0.25, 0.1 , 0.05])
And sure enough, amps2 is the sequence of amplitudes we started with.
This algorithm works, but it is slow. Solving a linear system of equations takes time proportional to \(n^3\), where \(n\) is the number of columns in \(M\). We can do better.
6.4. Orthogonal matrices#
One way to solve linear systems is by inverting matrices. The inverse of a matrix \(M\) is written \(M^{-1}\), and it has the property that \(M^{-1}M = I\). \(I\) is the identity matrix, which has the value 1 on all diagonal elements and 0 everywhere else.
So, to solve the equation \(y = Ma\), we can multiply both sides by \(M^{-1}\), which yields:
On the right side, we can replace \(M^{-1}M\) with \(I\):
If we multiply \(I\) by any vector \(a\), the result is \(a\), so:
This implies that if we can compute \(M^{-1}\) efficiently, we can find \(a\) with a simple matrix multiplication (using np.dot).
That takes time proportional to \(n^2\), which is better than \(n^3\).
Inverting a matrix is slow, in general, but some special cases are faster. In particular, if \(M\) is orthogonal, the inverse of \(M\) is just the transpose of \(M\), written \(M^T\). In NumPy transposing an array is a constant-time operation. It doesn’t actually move the elements of the array; instead, it creates a “view” that changes the way the elements are accessed.
Again, a matrix is orthogonal if its transpose is also its inverse; that is, \(M^T = M^{-1}\). That implies that \(M^TM = I\), which means we can check whether a matrix is orthogonal by computing \(M^TM\).
So let’s see what the matrix looks like in synthesize2.
In the previous example, \(M\) has 11,025 rows, so it might be a good idea to work with a smaller example:
# suppress scientific notation for small numbers
np.set_printoptions(precision=3, suppress=True)
amps = np.array([0.6, 0.25, 0.1, 0.05])
N = 4.0
time_unit = 0.001
ts = np.arange(N) / N * time_unit
max_freq = N / time_unit / 2
fs = np.arange(N) / N * max_freq
ys = synthesize2(amps, fs, ts)
amps is the same vector of amplitudes we saw before.
Since we have 4 frequency components, we’ll sample the signal at 4 points in time.
That way, \(M\) is square.
ts is a vector of equally spaced sample times in the range from 0 to 1 time unit.
I chose the time unit to be 1 millisecond, but it is an arbitrary choice, and we will see in a minute that it drops out of the computation anyway.
Since the frame rate is \(N\) samples per time unit, the Nyquist frequency is N / time_unit / 2, which is 2000 Hz in this example.
So fs is a vector of equally spaced frequencies between 0 and 2000 Hz.
With these values of ts and fs, the matrix, \(M\), is:
args = np.outer(ts, fs)
M = np.cos(PI2 * args)
M
array([[ 1. , 1. , 1. , 1. ],
[ 1. , 0.707, 0. , -0.707],
[ 1. , 0. , -1. , -0. ],
[ 1. , -0.707, -0. , 0.707]])
You might recognize 0.707 as an approximation of \(\sqrt{2}/2\), which is \(\cos \pi/4\). You also might notice that this matrix is symmetric, which means that the element at \((j, k)\) always equals the element at \((k, j)\). This implies that \(M\) is its own transpose; that is, \(M^T = M\).
But sadly, \(M\) is not orthogonal. If we compute \(M^TM\), we get:
M.T @ M
array([[ 4., 1., -0., 1.],
[ 1., 2., 1., -0.],
[-0., 1., 2., 1.],
[ 1., -0., 1., 2.]])
And that’s not the identity matrix.
6.5. DCT-IV#
But if we choose ts and fs carefully, we can make \(M\) orthogonal.
There are several ways to do it, which is why there are several versions of the discrete cosine transform (DCT).
One simple option is to shift ts and fs by a half unit.
This version is called DCT-IV, where “IV” is a roman numeral indicating that this is the fourth of eight versions of the DCT.
Here’s an updated version of test1:
N = 4.0
ts = (0.5 + np.arange(N)) / N
fs = (0.5 + np.arange(N)) / 2
If you compare this to the previous version, you’ll notice two changes.
First, I added 0.5 to ts and fs.
Second, I canceled out time_units, which simplifies the expression for fs.
With these values, \(M\) is
args = np.outer(ts, fs)
M = np.cos(PI2 * args)
M
array([[ 0.981, 0.831, 0.556, 0.195],
[ 0.831, -0.195, -0.981, -0.556],
[ 0.556, -0.981, 0.195, 0.831],
[ 0.195, -0.556, 0.831, -0.981]])
And \(M^TM\) is
M.T @ M
array([[ 2., -0., 0., 0.],
[-0., 2., -0., -0.],
[ 0., -0., 2., -0.],
[ 0., -0., -0., 2.]])
Some of the off-diagonal elements are displayed as -0, which means that the floating-point representation is a small negative number. This matrix is very close to \(2I\), which means \(M\) is almost orthogonal; it’s just off by a factor of 2. And for our purposes, that’s good enough.
Because \(M\) is symmetric and (almost) orthogonal, the inverse of \(M\) is just \(M/2\).
Now we can write a more efficient version of analyze:
def analyze2(ys, fs, ts):
args = np.outer(ts, fs)
M = np.cos(PI2 * args)
amps = np.dot(M, ys) / 2
return amps
Instead of using np.linalg.solve, we just multiply by \(M/2\).
Combining test2 and analyze2, we can write an implementation of DCT-IV:
def dct_iv(ys):
N = len(ys)
ts = (0.5 + np.arange(N)) / N
fs = (0.5 + np.arange(N)) / 2
args = np.outer(ts, fs)
M = np.cos(PI2 * args)
amps = np.dot(M, ys) / 2
return amps
Again, ys is the wave array.
We don’t have to pass ts and fs as parameters; dct_iv can figure them out based on N, the length of ys.
If we’ve got it right, this function should solve the analysis problem; that is, given ys it should be able to recover amps.
We can test it like this:
N = 4.0
ts = (0.5 + np.arange(N)) / N
fs = (0.5 + np.arange(N)) / 2
ys = synthesize2(amps, fs, ts)
amps2 = dct_iv(ys)
max(abs(amps - amps2))
np.float64(5.551115123125783e-17)
Starting with amps, we synthesize a wave array, then use dct_iv to compute amps2.
The biggest difference between amps and amps2 is about 1e-16, which is what we expect due to floating-point errors.
6.6. Inverse DCT#
Finally, notice that analyze2 and synthesize2 are almost identical.
The only difference is that analyze2 divides the result by 2. We can use this insight to compute the inverse DCT:
def inverse_dct_iv(amps):
return dct_iv(amps) * 2
inverse_dct_iv solves the synthesis problem: it takes the vector of amplitudes and returns the wave array, ys.
We can test it by starting with amps, applying inverse_dct_iv and dct_iv, and testing that we get back what we started with.
amps = [0.6, 0.25, 0.1, 0.05]
ys = inverse_dct_iv(amps)
amps2 = dct_iv(ys)
max(abs(amps - amps2))
np.float64(5.551115123125783e-17)
Again, the biggest difference is about 1e-16.
6.7. The Dct Class#
thinkdsp provides a Dct class that encapsulates the DCT in the same way the Spectrum class encapsulates the FFT.
To make a Dct object, you can invoke make_dct on a Wave.
from thinkdsp import TriangleSignal
signal = TriangleSignal(freq=400)
wave = signal.make_wave(duration=1.0, framerate=10000)
dct = wave.make_dct()
dct.plot()
decorate_freq()
The result is the DCT of a triangle wave at 400 Hz, as shown in the following figure. The values of the DCT can be positive or negative; a negative value in the DCT corresponds to a negated cosine or, equivalently, to a cosine shifted by 180 degrees.
make_dct uses DCT-II, which is the most common type of DCT, provided by scipy.fftpack.
import scipy.fftpack
# class Wave:
def make_dct(self):
N = len(self.ys)
hs = scipy.fftpack.dct(self.ys, type=2)
fs = (0.5 + np.arange(N)) / 2
return Dct(hs, fs, self.framerate)
The results from dct are stored in hs.
The corresponding frequencies, computed as in Section [dctiv]{reference-type=“ref” reference=“dctiv”}, are stored in fs.
And then both are used to initialize the Dct object.
Dct provides make_wave, which performs the inverse DCT.
We can test it like this:
wave2 = dct.make_wave()
max(abs(wave.ys-wave2.ys))
np.float64(7.771561172376096e-16)
The biggest difference between ys1 and ys2 is about 1e-16, which is what we expect due to floating-point errors.
make_wave uses scipy.fftpack.idct:
# class Dct
def make_wave(self):
n = len(self.hs)
ys = scipy.fftpack.idct(self.hs, type=2) / 2 / n
return Wave(ys, framerate=self.framerate)
By default, the inverse DCT doesn’t normalize the result, so we have to divide through by \(2N\).
6.8. Exercises#
6.9. Exercise 1#
In this chapter I claim that analyze1 takes time proportional to \(n^3\) and analyze2 takes time proportional to \(n^2\).
To see if that’s true, run them on a range of input sizes and time them.
In IPython, you can use the magic command %timeit.
If you plot run time versus input size on a log-log scale, you should get a straight line with slope 3 for analyze1 and slope 2 for analyze2.
You also might want to test dct_iv and scipy.fftpack.dct.
I’ll start with a noise signal and an array of power-of-two sizes
from thinkdsp import UncorrelatedGaussianNoise
signal = UncorrelatedGaussianNoise()
noise = signal.make_wave(duration=1.0, framerate=16384)
noise.ys.shape
(16384,)
The following function takes an array of results from a timing experiment, plots the results, and fits a straight line.
from scipy.stats import linregress
loglog = dict(xscale='log', yscale='log')
def plot_bests(ns, bests):
plt.plot(ns, bests)
decorate(**loglog)
x = np.log(ns)
y = np.log(bests)
t = linregress(x,y)
slope = t[0]
return slope
def analyze1(ys, fs, ts):
"""Analyze a mixture of cosines and return amplitudes.
Works for the general case where M is not orthogonal.
ys: wave array
fs: frequencies in Hz
ts: times where the signal was evaluated
returns: vector of amplitudes
"""
args = np.outer(ts, fs)
M = np.cos(PI2 * args)
amps = np.linalg.solve(M, ys)
return amps
def run_speed_test(ns, func):
results = []
for N in ns:
print(N)
ts = (0.5 + np.arange(N)) / N
freqs = (0.5 + np.arange(N)) / 2
ys = noise.ys[:N]
result = %timeit -r1 -o func(ys, freqs, ts)
results.append(result)
bests = [result.best for result in results]
return bests
Here are the results for analyze1.
ns = 2 ** np.arange(6, 13)
ns
array([ 64, 128, 256, 512, 1024, 2048, 4096])
bests = run_speed_test(ns, analyze1)
plot_bests(ns, bests)
64
212 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 10,000 loops each)
128
566 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 1,000 loops each)
256
2.59 ms ± 0 ns per loop (mean ± std. dev. of 1 run, 100 loops each)
512
22.5 ms ± 0 ns per loop (mean ± std. dev. of 1 run, 100 loops each)
1024
99.5 ms ± 0 ns per loop (mean ± std. dev. of 1 run, 10 loops each)
2048
671 ms ± 0 ns per loop (mean ± std. dev. of 1 run, 1 loop each)
4096
3.36 s ± 0 ns per loop (mean ± std. dev. of 1 run, 1 loop each)
np.float64(2.4122983175572688)
The estimated slope is close to 2, not 3, as expected.
One possibility is that the performance of np.linalg.solve is nearly quadratic in this range of array sizes.
Here are the results for analyze2:
def analyze2(ys, fs, ts):
"""Analyze a mixture of cosines and return amplitudes.
Assumes that fs and ts are chosen so that M is orthogonal.
ys: wave array
fs: frequencies in Hz
ts: times where the signal was evaluated
returns: vector of amplitudes
"""
args = np.outer(ts, fs)
M = np.cos(PI2 * args)
amps = np.dot(M, ys) / 2
return amps
bests2 = run_speed_test(ns, analyze2)
plot_bests(ns, bests2)
64
131 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 10,000 loops each)
128
344 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 1,000 loops each)
256
1.42 ms ± 0 ns per loop (mean ± std. dev. of 1 run, 1,000 loops each)
512
5.84 ms ± 0 ns per loop (mean ± std. dev. of 1 run, 100 loops each)
1024
45.8 ms ± 0 ns per loop (mean ± std. dev. of 1 run, 10 loops each)
2048
248 ms ± 0 ns per loop (mean ± std. dev. of 1 run, 1 loop each)
4096
912 ms ± 0 ns per loop (mean ± std. dev. of 1 run, 1 loop each)
np.float64(2.224443259340936)
The results for analyze2 fall in a straight line with the estimated slope close to 2, as expected.
Here are the results for the scipy.fftpack.dct
import scipy.fftpack
def scipy_dct(ys, freqs, ts):
return scipy.fftpack.dct(ys, type=3)
bests3 = run_speed_test(ns, scipy_dct)
plot_bests(ns, bests3)
64
7.55 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 100,000 loops each)
128
8 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 100,000 loops each)
256
8.96 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 100,000 loops each)
512
11.1 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 100,000 loops each)
1024
14.9 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 100,000 loops each)
2048
24 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 10,000 loops each)
4096
41.2 μs ± 0 ns per loop (mean ± std. dev. of 1 run, 10,000 loops each)
np.float64(0.4016271288214343)
This implementation of dct is even faster. The line is curved, which means either we haven’t seen the asymptotic behavior yet, or the asymptotic behavior is not a simple exponent of \(n\). In fact, as we’ll see soon, the run time is proportional to \(n \log n\).
The following figure shows all three curves on the same axes.
plt.plot(ns, bests, label='analyze1')
plt.plot(ns, bests2, label='analyze2')
plt.plot(ns, bests3, label='fftpack.dct')
decorate(xlabel='Wave length (N)', ylabel='Time (s)', **loglog)
6.10. Exercise 2#
One of the major applications of the DCT is compression for both sound and images. In its simplest form, DCT-based compression works like this:
Break a long signal into segments.
Compute the DCT of each segment.
Identify frequency components with amplitudes so low they are inaudible, and remove them. Store only the frequencies and amplitudes that remain.
To play back the signal, load the frequencies and amplitudes for each segment and apply the inverse DCT.
Implement a version of this algorithm and apply it to a recording of music or speech. How many components can you eliminate before the difference is perceptible?
thinkdsp provides a class, Dct that is similar to a Spectrum, but which uses DCT instead of FFT.
As an example, I’ll use a recording of a saxophone:
download('https://github.com/AllenDowney/ThinkDSP/raw/master/code/100475__iluppai__saxophone-weep.wav')
Downloaded 100475__iluppai__saxophone-weep.wav
from thinkdsp import read_wave
wave = read_wave('100475__iluppai__saxophone-weep.wav')
wave.make_audio()
Here’s a short segment:
segment = wave.segment(start=1.2, duration=0.5)
segment.normalize()
segment.make_audio()
And here’s the DCT of that segment:
seg_dct = segment.make_dct()
seg_dct.plot(high=4000)
decorate_freq()
There are only a few harmonics with substantial amplitude, and many entries near zero.
The following function takes a DCT and sets elements below thresh to 0.
def compress(dct, thresh=1):
count = 0
for i, amp in enumerate(dct.amps):
if np.abs(amp) < thresh:
dct.hs[i] = 0
count += 1
n = len(dct.amps)
print(count, n, 100 * count / n, sep='\t')
If we apply it to the segment, we can eliminate more than 90% of the elements:
seg_dct = segment.make_dct()
compress(seg_dct, thresh=10)
seg_dct.plot(high=4000)
20457 22050 92.77551020408163
And the result sounds the same (at least to me):
seg2 = seg_dct.make_wave()
seg2.make_audio()
To compress a longer segment, we can make a DCT spectrogram.
The following function is similar to wave.make_spectrogram except that it uses the DCT.
from thinkdsp import Spectrogram
def make_dct_spectrogram(wave, seg_length):
"""Computes the DCT spectrogram of the wave.
seg_length: number of samples in each segment
returns: Spectrogram
"""
window = np.hamming(seg_length)
i, j = 0, seg_length
step = seg_length // 2
# map from time to Spectrum
spec_map = {}
while j < len(wave.ys):
segment = wave.slice(i, j)
segment.window(window)
# the nominal time for this segment is the midpoint
t = (segment.start + segment.end) / 2
spec_map[t] = segment.make_dct()
i += step
j += step
return Spectrogram(spec_map, seg_length)
Now we can make a DCT spectrogram and apply compress to each segment:
spectro = make_dct_spectrogram(wave, seg_length=1024)
for t, dct in sorted(spectro.spec_map.items()):
compress(dct, thresh=0.2)
1018 1024 99.4140625
1016 1024 99.21875
1014 1024 99.0234375
1017 1024 99.31640625
1016 1024 99.21875
1017 1024 99.31640625
1016 1024 99.21875
1020 1024 99.609375
1014 1024 99.0234375
1005 1024 98.14453125
1009 1024 98.53515625
1015 1024 99.12109375
1015 1024 99.12109375
1016 1024 99.21875
1016 1024 99.21875
1015 1024 99.12109375
1017 1024 99.31640625
1020 1024 99.609375
1013 1024 98.92578125
1017 1024 99.31640625
1013 1024 98.92578125
1017 1024 99.31640625
1018 1024 99.4140625
1015 1024 99.12109375
1013 1024 98.92578125
794 1024 77.5390625
785 1024 76.66015625
955 1024 93.26171875
995 1024 97.16796875
992 1024 96.875
976 1024 95.3125
925 1024 90.33203125
802 1024 78.3203125
836 1024 81.640625
850 1024 83.0078125
882 1024 86.1328125
883 1024 86.23046875
891 1024 87.01171875
901 1024 87.98828125
902 1024 88.0859375
900 1024 87.890625
900 1024 87.890625
894 1024 87.3046875
904 1024 88.28125
901 1024 87.98828125
915 1024 89.35546875
913 1024 89.16015625
899 1024 87.79296875
905 1024 88.37890625
905 1024 88.37890625
888 1024 86.71875
898 1024 87.6953125
879 1024 85.83984375
893 1024 87.20703125
893 1024 87.20703125
882 1024 86.1328125
874 1024 85.3515625
876 1024 85.546875
864 1024 84.375
879 1024 85.83984375
869 1024 84.86328125
872 1024 85.15625
871 1024 85.05859375
878 1024 85.7421875
872 1024 85.15625
859 1024 83.88671875
879 1024 85.83984375
889 1024 86.81640625
872 1024 85.15625
837 1024 81.73828125
842 1024 82.2265625
825 1024 80.56640625
839 1024 81.93359375
796 1024 77.734375
792 1024 77.34375
769 1024 75.09765625
836 1024 81.640625
919 1024 89.74609375
913 1024 89.16015625
942 1024 91.9921875
837 1024 81.73828125
739 1024 72.16796875
737 1024 71.97265625
726 1024 70.8984375
728 1024 71.09375
733 1024 71.58203125
717 1024 70.01953125
716 1024 69.921875
676 1024 66.015625
712 1024 69.53125
697 1024 68.06640625
718 1024 70.1171875
717 1024 70.01953125
718 1024 70.1171875
681 1024 66.50390625
707 1024 69.04296875
691 1024 67.48046875
681 1024 66.50390625
709 1024 69.23828125
684 1024 66.796875
743 1024 72.55859375
710 1024 69.3359375
712 1024 69.53125
714 1024 69.7265625
719 1024 70.21484375
708 1024 69.140625
725 1024 70.80078125
700 1024 68.359375
726 1024 70.8984375
716 1024 69.921875
725 1024 70.80078125
692 1024 67.578125
675 1024 65.91796875
747 1024 72.94921875
741 1024 72.36328125
730 1024 71.2890625
701 1024 68.45703125
721 1024 70.41015625
747 1024 72.94921875
725 1024 70.80078125
744 1024 72.65625
720 1024 70.3125
716 1024 69.921875
723 1024 70.60546875
721 1024 70.41015625
734 1024 71.6796875
730 1024 71.2890625
718 1024 70.1171875
730 1024 71.2890625
723 1024 70.60546875
749 1024 73.14453125
727 1024 70.99609375
728 1024 71.09375
746 1024 72.8515625
739 1024 72.16796875
757 1024 73.92578125
741 1024 72.36328125
751 1024 73.33984375
775 1024 75.68359375
749 1024 73.14453125
768 1024 75.0
763 1024 74.51171875
771 1024 75.29296875
758 1024 74.0234375
745 1024 72.75390625
756 1024 73.828125
744 1024 72.65625
743 1024 72.55859375
757 1024 73.92578125
779 1024 76.07421875
760 1024 74.21875
770 1024 75.1953125
759 1024 74.12109375
737 1024 71.97265625
739 1024 72.16796875
751 1024 73.33984375
762 1024 74.4140625
754 1024 73.6328125
811 1024 79.19921875
899 1024 87.79296875
832 1024 81.25
800 1024 78.125
756 1024 73.828125
748 1024 73.046875
727 1024 70.99609375
744 1024 72.65625
725 1024 70.80078125
720 1024 70.3125
755 1024 73.73046875
737 1024 71.97265625
766 1024 74.8046875
747 1024 72.94921875
743 1024 72.55859375
727 1024 70.99609375
726 1024 70.8984375
746 1024 72.8515625
764 1024 74.609375
751 1024 73.33984375
734 1024 71.6796875
741 1024 72.36328125
760 1024 74.21875
750 1024 73.2421875
784 1024 76.5625
730 1024 71.2890625
757 1024 73.92578125
761 1024 74.31640625
734 1024 71.6796875
744 1024 72.65625
757 1024 73.92578125
714 1024 69.7265625
740 1024 72.265625
738 1024 72.0703125
763 1024 74.51171875
766 1024 74.8046875
745 1024 72.75390625
751 1024 73.33984375
759 1024 74.12109375
756 1024 73.828125
756 1024 73.828125
756 1024 73.828125
755 1024 73.73046875
746 1024 72.8515625
756 1024 73.828125
738 1024 72.0703125
757 1024 73.92578125
764 1024 74.609375
765 1024 74.70703125
762 1024 74.4140625
768 1024 75.0
773 1024 75.48828125
782 1024 76.3671875
773 1024 75.48828125
766 1024 74.8046875
755 1024 73.73046875
766 1024 74.8046875
772 1024 75.390625
810 1024 79.1015625
739 1024 72.16796875
717 1024 70.01953125
722 1024 70.5078125
739 1024 72.16796875
725 1024 70.80078125
736 1024 71.875
759 1024 74.12109375
769 1024 75.09765625
749 1024 73.14453125
710 1024 69.3359375
748 1024 73.046875
720 1024 70.3125
732 1024 71.484375
721 1024 70.41015625
734 1024 71.6796875
763 1024 74.51171875
747 1024 72.94921875
754 1024 73.6328125
755 1024 73.73046875
764 1024 74.609375
801 1024 78.22265625
768 1024 75.0
780 1024 76.171875
773 1024 75.48828125
764 1024 74.609375
775 1024 75.68359375
740 1024 72.265625
794 1024 77.5390625
796 1024 77.734375
769 1024 75.09765625
751 1024 73.33984375
782 1024 76.3671875
758 1024 74.0234375
777 1024 75.87890625
794 1024 77.5390625
784 1024 76.5625
788 1024 76.953125
773 1024 75.48828125
783 1024 76.46484375
784 1024 76.5625
785 1024 76.66015625
806 1024 78.7109375
807 1024 78.80859375
797 1024 77.83203125
785 1024 76.66015625
794 1024 77.5390625
766 1024 74.8046875
790 1024 77.1484375
746 1024 72.8515625
762 1024 74.4140625
813 1024 79.39453125
801 1024 78.22265625
782 1024 76.3671875
776 1024 75.78125
755 1024 73.73046875
780 1024 76.171875
784 1024 76.5625
805 1024 78.61328125
791 1024 77.24609375
803 1024 78.41796875
799 1024 78.02734375
795 1024 77.63671875
797 1024 77.83203125
806 1024 78.7109375
781 1024 76.26953125
795 1024 77.63671875
797 1024 77.83203125
893 1024 87.20703125
775 1024 75.68359375
787 1024 76.85546875
746 1024 72.8515625
767 1024 74.90234375
749 1024 73.14453125
749 1024 73.14453125
738 1024 72.0703125
736 1024 71.875
747 1024 72.94921875
760 1024 74.21875
737 1024 71.97265625
752 1024 73.4375
756 1024 73.828125
772 1024 75.390625
740 1024 72.265625
737 1024 71.97265625
766 1024 74.8046875
791 1024 77.24609375
765 1024 74.70703125
771 1024 75.29296875
786 1024 76.7578125
770 1024 75.1953125
761 1024 74.31640625
765 1024 74.70703125
756 1024 73.828125
758 1024 74.0234375
765 1024 74.70703125
785 1024 76.66015625
769 1024 75.09765625
781 1024 76.26953125
792 1024 77.34375
798 1024 77.9296875
809 1024 79.00390625
778 1024 75.9765625
782 1024 76.3671875
776 1024 75.78125
791 1024 77.24609375
794 1024 77.5390625
783 1024 76.46484375
771 1024 75.29296875
792 1024 77.34375
785 1024 76.66015625
812 1024 79.296875
809 1024 79.00390625
799 1024 78.02734375
798 1024 77.9296875
803 1024 78.41796875
800 1024 78.125
805 1024 78.61328125
803 1024 78.41796875
799 1024 78.02734375
802 1024 78.3203125
804 1024 78.515625
809 1024 79.00390625
784 1024 76.5625
791 1024 77.24609375
814 1024 79.4921875
788 1024 76.953125
816 1024 79.6875
810 1024 79.1015625
820 1024 80.078125
823 1024 80.37109375
813 1024 79.39453125
799 1024 78.02734375
807 1024 78.80859375
799 1024 78.02734375
789 1024 77.05078125
813 1024 79.39453125
819 1024 79.98046875
809 1024 79.00390625
784 1024 76.5625
809 1024 79.00390625
810 1024 79.1015625
785 1024 76.66015625
838 1024 81.8359375
821 1024 80.17578125
822 1024 80.2734375
800 1024 78.125
815 1024 79.58984375
827 1024 80.76171875
820 1024 80.078125
792 1024 77.34375
818 1024 79.8828125
813 1024 79.39453125
824 1024 80.46875
795 1024 77.63671875
788 1024 76.953125
796 1024 77.734375
802 1024 78.3203125
800 1024 78.125
796 1024 77.734375
823 1024 80.37109375
804 1024 78.515625
811 1024 79.19921875
808 1024 78.90625
815 1024 79.58984375
812 1024 79.296875
822 1024 80.2734375
793 1024 77.44140625
803 1024 78.41796875
806 1024 78.7109375
812 1024 79.296875
796 1024 77.734375
804 1024 78.515625
807 1024 78.80859375
821 1024 80.17578125
793 1024 77.44140625
799 1024 78.02734375
810 1024 79.1015625
818 1024 79.8828125
813 1024 79.39453125
825 1024 80.56640625
804 1024 78.515625
821 1024 80.17578125
809 1024 79.00390625
828 1024 80.859375
813 1024 79.39453125
838 1024 81.8359375
836 1024 81.640625
818 1024 79.8828125
808 1024 78.90625
819 1024 79.98046875
820 1024 80.078125
814 1024 79.4921875
901 1024 87.98828125
894 1024 87.3046875
888 1024 86.71875
780 1024 76.171875
773 1024 75.48828125
750 1024 73.2421875
750 1024 73.2421875
730 1024 71.2890625
761 1024 74.31640625
775 1024 75.68359375
782 1024 76.3671875
788 1024 76.953125
748 1024 73.046875
752 1024 73.4375
771 1024 75.29296875
746 1024 72.8515625
778 1024 75.9765625
777 1024 75.87890625
760 1024 74.21875
734 1024 71.6796875
711 1024 69.43359375
754 1024 73.6328125
745 1024 72.75390625
758 1024 74.0234375
744 1024 72.65625
755 1024 73.73046875
749 1024 73.14453125
723 1024 70.60546875
784 1024 76.5625
761 1024 74.31640625
758 1024 74.0234375
709 1024 69.23828125
769 1024 75.09765625
773 1024 75.48828125
769 1024 75.09765625
756 1024 73.828125
747 1024 72.94921875
787 1024 76.85546875
770 1024 75.1953125
749 1024 73.14453125
769 1024 75.09765625
748 1024 73.046875
761 1024 74.31640625
759 1024 74.12109375
775 1024 75.68359375
756 1024 73.828125
774 1024 75.5859375
776 1024 75.78125
760 1024 74.21875
783 1024 76.46484375
744 1024 72.65625
766 1024 74.8046875
761 1024 74.31640625
788 1024 76.953125
774 1024 75.5859375
753 1024 73.53515625
754 1024 73.6328125
765 1024 74.70703125
736 1024 71.875
782 1024 76.3671875
768 1024 75.0
778 1024 75.9765625
767 1024 74.90234375
774 1024 75.5859375
772 1024 75.390625
769 1024 75.09765625
774 1024 75.5859375
776 1024 75.78125
796 1024 77.734375
762 1024 74.4140625
766 1024 74.8046875
765 1024 74.70703125
783 1024 76.46484375
770 1024 75.1953125
799 1024 78.02734375
779 1024 76.07421875
774 1024 75.5859375
791 1024 77.24609375
797 1024 77.83203125
781 1024 76.26953125
754 1024 73.6328125
790 1024 77.1484375
790 1024 77.1484375
801 1024 78.22265625
783 1024 76.46484375
787 1024 76.85546875
805 1024 78.61328125
758 1024 74.0234375
785 1024 76.66015625
788 1024 76.953125
806 1024 78.7109375
818 1024 79.8828125
776 1024 75.78125
807 1024 78.80859375
802 1024 78.3203125
782 1024 76.3671875
812 1024 79.296875
803 1024 78.41796875
803 1024 78.41796875
787 1024 76.85546875
799 1024 78.02734375
786 1024 76.7578125
813 1024 79.39453125
813 1024 79.39453125
813 1024 79.39453125
803 1024 78.41796875
815 1024 79.58984375
792 1024 77.34375
807 1024 78.80859375
829 1024 80.95703125
797 1024 77.83203125
814 1024 79.4921875
793 1024 77.44140625
802 1024 78.3203125
775 1024 75.68359375
816 1024 79.6875
804 1024 78.515625
808 1024 78.90625
809 1024 79.00390625
814 1024 79.4921875
808 1024 78.90625
823 1024 80.37109375
811 1024 79.19921875
806 1024 78.7109375
819 1024 79.98046875
805 1024 78.61328125
826 1024 80.6640625
826 1024 80.6640625
807 1024 78.80859375
818 1024 79.8828125
818 1024 79.8828125
812 1024 79.296875
816 1024 79.6875
815 1024 79.58984375
827 1024 80.76171875
830 1024 81.0546875
852 1024 83.203125
827 1024 80.76171875
834 1024 81.4453125
835 1024 81.54296875
835 1024 81.54296875
829 1024 80.95703125
822 1024 80.2734375
818 1024 79.8828125
827 1024 80.76171875
834 1024 81.4453125
829 1024 80.95703125
846 1024 82.6171875
829 1024 80.95703125
829 1024 80.95703125
833 1024 81.34765625
837 1024 81.73828125
837 1024 81.73828125
815 1024 79.58984375
834 1024 81.4453125
833 1024 81.34765625
840 1024 82.03125
855 1024 83.49609375
853 1024 83.30078125
853 1024 83.30078125
846 1024 82.6171875
852 1024 83.203125
856 1024 83.59375
859 1024 83.88671875
851 1024 83.10546875
845 1024 82.51953125
874 1024 85.3515625
861 1024 84.08203125
877 1024 85.64453125
853 1024 83.30078125
861 1024 84.08203125
859 1024 83.88671875
866 1024 84.5703125
868 1024 84.765625
870 1024 84.9609375
856 1024 83.59375
859 1024 83.88671875
864 1024 84.375
864 1024 84.375
876 1024 85.546875
872 1024 85.15625
872 1024 85.15625
863 1024 84.27734375
859 1024 83.88671875
878 1024 85.7421875
860 1024 83.984375
864 1024 84.375
875 1024 85.44921875
862 1024 84.1796875
867 1024 84.66796875
867 1024 84.66796875
864 1024 84.375
864 1024 84.375
876 1024 85.546875
875 1024 85.44921875
860 1024 83.984375
865 1024 84.47265625
881 1024 86.03515625
867 1024 84.66796875
869 1024 84.86328125
873 1024 85.25390625
869 1024 84.86328125
873 1024 85.25390625
873 1024 85.25390625
862 1024 84.1796875
865 1024 84.47265625
871 1024 85.05859375
869 1024 84.86328125
871 1024 85.05859375
866 1024 84.5703125
877 1024 85.64453125
861 1024 84.08203125
881 1024 86.03515625
882 1024 86.1328125
874 1024 85.3515625
875 1024 85.44921875
866 1024 84.5703125
870 1024 84.9609375
883 1024 86.23046875
870 1024 84.9609375
871 1024 85.05859375
877 1024 85.64453125
866 1024 84.5703125
877 1024 85.64453125
863 1024 84.27734375
873 1024 85.25390625
871 1024 85.05859375
883 1024 86.23046875
862 1024 84.1796875
853 1024 83.30078125
858 1024 83.7890625
857 1024 83.69140625
855 1024 83.49609375
847 1024 82.71484375
837 1024 81.73828125
850 1024 83.0078125
864 1024 84.375
879 1024 85.83984375
883 1024 86.23046875
871 1024 85.05859375
888 1024 86.71875
881 1024 86.03515625
830 1024 81.0546875
870 1024 84.9609375
877 1024 85.64453125
886 1024 86.5234375
863 1024 84.27734375
871 1024 85.05859375
886 1024 86.5234375
871 1024 85.05859375
896 1024 87.5
872 1024 85.15625
870 1024 84.9609375
877 1024 85.64453125
863 1024 84.27734375
886 1024 86.5234375
898 1024 87.6953125
884 1024 86.328125
908 1024 88.671875
878 1024 85.7421875
865 1024 84.47265625
864 1024 84.375
888 1024 86.71875
870 1024 84.9609375
862 1024 84.1796875
866 1024 84.5703125
889 1024 86.81640625
879 1024 85.83984375
884 1024 86.328125
880 1024 85.9375
876 1024 85.546875
864 1024 84.375
877 1024 85.64453125
858 1024 83.7890625
894 1024 87.3046875
890 1024 86.9140625
893 1024 87.20703125
891 1024 87.01171875
896 1024 87.5
892 1024 87.109375
906 1024 88.4765625
878 1024 85.7421875
893 1024 87.20703125
898 1024 87.6953125
888 1024 86.71875
903 1024 88.18359375
911 1024 88.96484375
911 1024 88.96484375
901 1024 87.98828125
909 1024 88.76953125
911 1024 88.96484375
921 1024 89.94140625
922 1024 90.0390625
916 1024 89.453125
923 1024 90.13671875
928 1024 90.625
920 1024 89.84375
922 1024 90.0390625
915 1024 89.35546875
930 1024 90.8203125
914 1024 89.2578125
917 1024 89.55078125
918 1024 89.6484375
921 1024 89.94140625
921 1024 89.94140625
937 1024 91.50390625
931 1024 90.91796875
923 1024 90.13671875
921 1024 89.94140625
934 1024 91.2109375
930 1024 90.8203125
933 1024 91.11328125
932 1024 91.015625
930 1024 90.8203125
930 1024 90.8203125
933 1024 91.11328125
933 1024 91.11328125
949 1024 92.67578125
941 1024 91.89453125
945 1024 92.28515625
936 1024 91.40625
956 1024 93.359375
948 1024 92.578125
936 1024 91.40625
941 1024 91.89453125
949 1024 92.67578125
941 1024 91.89453125
940 1024 91.796875
951 1024 92.87109375
941 1024 91.89453125
941 1024 91.89453125
930 1024 90.8203125
930 1024 90.8203125
924 1024 90.234375
919 1024 89.74609375
911 1024 88.96484375
934 1024 91.2109375
892 1024 87.109375
929 1024 90.72265625
922 1024 90.0390625
927 1024 90.52734375
917 1024 89.55078125
856 1024 83.59375
835 1024 81.54296875
852 1024 83.203125
870 1024 84.9609375
878 1024 85.7421875
872 1024 85.15625
894 1024 87.3046875
865 1024 84.47265625
889 1024 86.81640625
871 1024 85.05859375
873 1024 85.25390625
864 1024 84.375
859 1024 83.88671875
867 1024 84.66796875
833 1024 81.34765625
853 1024 83.30078125
874 1024 85.3515625
843 1024 82.32421875
848 1024 82.8125
844 1024 82.421875
817 1024 79.78515625
865 1024 84.47265625
807 1024 78.80859375
752 1024 73.4375
775 1024 75.68359375
772 1024 75.390625
778 1024 75.9765625
767 1024 74.90234375
784 1024 76.5625
800 1024 78.125
807 1024 78.80859375
826 1024 80.6640625
805 1024 78.61328125
788 1024 76.953125
820 1024 80.078125
809 1024 79.00390625
803 1024 78.41796875
799 1024 78.02734375
806 1024 78.7109375
839 1024 81.93359375
846 1024 82.6171875
914 1024 89.2578125
888 1024 86.71875
839 1024 81.93359375
836 1024 81.640625
830 1024 81.0546875
845 1024 82.51953125
828 1024 80.859375
834 1024 81.4453125
854 1024 83.3984375
847 1024 82.71484375
846 1024 82.6171875
845 1024 82.51953125
863 1024 84.27734375
867 1024 84.66796875
855 1024 83.49609375
844 1024 82.421875
864 1024 84.375
865 1024 84.47265625
860 1024 83.984375
868 1024 84.765625
871 1024 85.05859375
868 1024 84.765625
857 1024 83.69140625
885 1024 86.42578125
908 1024 88.671875
872 1024 85.15625
848 1024 82.8125
813 1024 79.39453125
763 1024 74.51171875
761 1024 74.31640625
779 1024 76.07421875
783 1024 76.46484375
776 1024 75.78125
784 1024 76.5625
811 1024 79.19921875
812 1024 79.296875
777 1024 75.87890625
794 1024 77.5390625
794 1024 77.5390625
813 1024 79.39453125
805 1024 78.61328125
828 1024 80.859375
796 1024 77.734375
806 1024 78.7109375
817 1024 79.78515625
840 1024 82.03125
786 1024 76.7578125
818 1024 79.8828125
832 1024 81.25
835 1024 81.54296875
768 1024 75.0
841 1024 82.12890625
833 1024 81.34765625
836 1024 81.640625
824 1024 80.46875
830 1024 81.0546875
837 1024 81.73828125
837 1024 81.73828125
858 1024 83.7890625
847 1024 82.71484375
870 1024 84.9609375
866 1024 84.5703125
843 1024 82.32421875
867 1024 84.66796875
849 1024 82.91015625
869 1024 84.86328125
860 1024 83.984375
862 1024 84.1796875
846 1024 82.6171875
854 1024 83.3984375
871 1024 85.05859375
861 1024 84.08203125
882 1024 86.1328125
892 1024 87.109375
877 1024 85.64453125
895 1024 87.40234375
887 1024 86.62109375
873 1024 85.25390625
894 1024 87.3046875
890 1024 86.9140625
878 1024 85.7421875
894 1024 87.3046875
879 1024 85.83984375
890 1024 86.9140625
895 1024 87.40234375
889 1024 86.81640625
896 1024 87.5
898 1024 87.6953125
901 1024 87.98828125
879 1024 85.83984375
890 1024 86.9140625
888 1024 86.71875
917 1024 89.55078125
902 1024 88.0859375
921 1024 89.94140625
915 1024 89.35546875
927 1024 90.52734375
923 1024 90.13671875
928 1024 90.625
923 1024 90.13671875
914 1024 89.2578125
918 1024 89.6484375
927 1024 90.52734375
926 1024 90.4296875
919 1024 89.74609375
916 1024 89.453125
928 1024 90.625
916 1024 89.453125
933 1024 91.11328125
925 1024 90.33203125
930 1024 90.8203125
930 1024 90.8203125
934 1024 91.2109375
933 1024 91.11328125
935 1024 91.30859375
939 1024 91.69921875
934 1024 91.2109375
938 1024 91.6015625
944 1024 92.1875
937 1024 91.50390625
937 1024 91.50390625
935 1024 91.30859375
937 1024 91.50390625
937 1024 91.50390625
954 1024 93.1640625
940 1024 91.796875
942 1024 91.9921875
955 1024 93.26171875
949 1024 92.67578125
941 1024 91.89453125
947 1024 92.48046875
940 1024 91.796875
943 1024 92.08984375
946 1024 92.3828125
962 1024 93.9453125
954 1024 93.1640625
956 1024 93.359375
957 1024 93.45703125
962 1024 93.9453125
960 1024 93.75
944 1024 92.1875
969 1024 94.62890625
969 1024 94.62890625
969 1024 94.62890625
968 1024 94.53125
970 1024 94.7265625
969 1024 94.62890625
975 1024 95.21484375
963 1024 94.04296875
965 1024 94.23828125
975 1024 95.21484375
969 1024 94.62890625
966 1024 94.3359375
969 1024 94.62890625
984 1024 96.09375
981 1024 95.80078125
977 1024 95.41015625
982 1024 95.8984375
980 1024 95.703125
980 1024 95.703125
981 1024 95.80078125
982 1024 95.8984375
981 1024 95.80078125
987 1024 96.38671875
982 1024 95.8984375
982 1024 95.8984375
977 1024 95.41015625
982 1024 95.8984375
986 1024 96.2890625
984 1024 96.09375
986 1024 96.2890625
987 1024 96.38671875
996 1024 97.265625
995 1024 97.16796875
997 1024 97.36328125
999 1024 97.55859375
1002 1024 97.8515625
999 1024 97.55859375
1004 1024 98.046875
1002 1024 97.8515625
999 1024 97.55859375
1000 1024 97.65625
1007 1024 98.33984375
1010 1024 98.6328125
1012 1024 98.828125
1004 1024 98.046875
1000 1024 97.65625
1008 1024 98.4375
1005 1024 98.14453125
1011 1024 98.73046875
1011 1024 98.73046875
1009 1024 98.53515625
1005 1024 98.14453125
1007 1024 98.33984375
1010 1024 98.6328125
1010 1024 98.6328125
1005 1024 98.14453125
1008 1024 98.4375
1012 1024 98.828125
1009 1024 98.53515625
1012 1024 98.828125
1013 1024 98.92578125
1010 1024 98.6328125
1012 1024 98.828125
1014 1024 99.0234375
1016 1024 99.21875
1010 1024 98.6328125
1014 1024 99.0234375
1015 1024 99.12109375
1012 1024 98.828125
1019 1024 99.51171875
1015 1024 99.12109375
1016 1024 99.21875
1019 1024 99.51171875
1019 1024 99.51171875
1016 1024 99.21875
1018 1024 99.4140625
1018 1024 99.4140625
In most segments, the compression is 75-85%.
To hear what it sounds like, we can convert the spectrogram back to a wave and play it.
wave2 = spectro.make_wave()
wave2.make_audio()
And here’s the original again for comparison.
wave.make_audio()
As an experiment, you might try increasing thresh to see when the effect of compression becomes audible (to you).
Also, you might try compressing a signal with some noisy elements, like cymbals.
Think DSP: Digital Signal Processing in Python, 2rd Edition
Copyright 2024 Allen B. Downey
License: Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International